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Numerical  simulation  of  transport  phenomena  in  a  six-cell  PEM  fuel  cell  stack  is  performed  to  adjust  the 
gas  channel  size  of  individual  cells  to  obtain  evenly  distributed  cell  voltages.  Commercial  software,  CFD- 
ACE+,  is  used  as  the  solution  code,  and  the  effects  of  various  combinations  of  geometric  parameters  of 
the  channels  are  investigated,  including  the  channel  width  ratios  (A)  and  the  heights  of  the  gas  channels 
(ft).  Furthermore,  a  six-cell  fuel  cell  stack  was  assembled  for  experiments,  and  the  polarization  curves 
of  the  stack  were  measured.  With  the  help  of  the  experimental  data,  the  values  of  the  physical  and 
electrochemical  parameters  adopted  in  the  computation  model  are  determined,  and  the  computation 
model  is  then  used  in  numerical  simulation.  Under  the  assumption  of  no  dry  out  in  the  cells,  an  even 
voltage  distribution  can  be  obtained  if  the  distribution  of  reactant  gas  is  uniform  among  the  cells.  In  this 
study,  better  combinations  of  the  channel  width  ratios  and  the  heights  of  gas  channels  for  individual 
cells  are  proposed  to  yield  a  uniform  distribution  of  the  reactant  gases  in  the  cells.  An  adjustment  of  the 
geometric  parameters  of  gas  channels  that  leads  to  evenly  distributed  voltages  is  attempted,  and  a  16.5% 
increase  in  the  performance  of  the  fuel  cell  stack  is  observed  simply  by  adjusting  the  sizes  of  the  gas 
channels. 

©  2009  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

A  single-cell  proton  exchange  membrane  fuel  cell  (PEMFC)  is 
equipped  with  a  polymer  electrolyte  membrane,  which  is  sand¬ 
wiched  between  two  electrodes  with  two  gas  diffusion  layers  to 
form  a  membrane-electrode  assembly  (MEA).  In  a  typical  fuel  cell 
structure  [1],  the  MEA  is  placed  between  two  carbon  plates  with 
machined  groves  on  the  inner  faces  that  are  used  as  flow  chan¬ 
nels  to  provide  sufficient  amount  of  fuel  and  oxidant.  Hydrogen 
gas  is  supplied  as  a  fuel  to  the  flow  channels  facing  the  gas  dif¬ 
fusion  layer  at  the  anode  side,  and  oxygen  gas,  or  air,  is  supplied 
to  the  flow  channels  facing  the  gas  diffusion  layer  at  the  cathode 
side.  Both  the  anode  and  the  cathode  sides  contain  a  catalyst  to 
speed  up  the  electro-chemical  reaction  processes.  The  fuel  is  oxi¬ 
dized  on  the  anode  side,  and  free  electrons  are  conducted  through 
the  load  in  an  external  circuit  to  the  cathode  side.  The  polymer 
electrolyte  membrane  prevents  fuel  and  oxidant  molecules  from 
mixing,  and,  therefore,  it  prevents  direct  combustion.  However,  the 
ions  are  allowed  to  migrate  through  the  polymer  electrolyte  mem¬ 
brane  toward  the  cathode  side  where  the  oxidant  and  the  electrons 
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are  provided.  The  electro-chemical  reaction  at  the  cathode  side 
produces  water  and  completes  the  circuit. 

A  number  of  theoretical  studies  have  considered  modeling  the 
electrochemical  reaction  and  transport  phenomena  in  a  single-cell 
PEMFC.  For  example,  a  three-dimensional  numerical  analysis  of 
a  straight  channel  in  a  PEM  fuel  cell  was  performed  by  Dutta  et 
al.  [2]  using  commercial  computational  fluid  dynamics  software. 
The  authors  found  that  membrane  thickness,  cell  voltage  and  cur¬ 
rent  density  can  affect  water  transport  across  the  membrane.  Hu 
et  al.  [3]  reported  a  numerical  signal  phase  model  that  includes 
a  three-dimensional  solution  to  the  fluid  flow,  species  transport 
and  electrochemical  reaction  in  PEM  fuel  cells  with  a  conventional 
and  an  interdigitated  flow  field.  They  found  that  the  interdigitated 
flow  field  has  better  performance  because  of  its  advantage  in  mass- 
transport  ability;  however,  it  has  a  larger  pressure  loss  for  its  flow 
through  the  electrode.  You  and  Liu  [4]  considered  liquid  water  accu¬ 
mulation  in  the  straight  channels  on  the  cathode  side  and  used 
a  multi-phase  model  to  obtain  the  simulation  results.  Recently, 
Lin  and  coworkers  [5,6]  investigated  optimization  of  the  proton 
exchange  membrane  fuel  cell  based  on  a  three-dimensional  numer¬ 
ical  model.  Optimization  methods  were  introduced  into  the  model 
to  determine  the  optimal  combination  of  the  design  parameters, 
including  the  channel  width  ratio,  the  porosity  of  gas  diffusion  layer 
(GDL),  and  the  porosity  of  the  catalyst  layer. 
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Nomenclature 

a  water  activity 

D  mass  diffusivity  (m2  s-1 ) 

F  Faraday  constant  (96,487  C  mol-1 ) 

h  height  of  gas  channel  (m) 

i  local  electric  current  density  (A  m-3 ) 

i  local  ionic  current  density  (A  m-3 ) 

I  electric  current  density  (A  m-2 ) 

j  transfer  current  density  (A  m-3 ) 

k  permeability  (m2) 

lc  gas  channel  width  (m) 

lL  gas  channel  length  (m) 

Ip  width  of  a  module  (m) 

M  molecular  weight  (kg  kmol-1 ) 

n  coordinate  in  the  direction  normal  to  the  surface 

P  pressure  (Pa) 

P  power  density  ( W  m-2 ) 

r  concentration  parameter 

R  universal  gas  constant  (8.314J  mol-1  K-1 ) 

5  source  term  in  species  equation 

T  temperature  (K) 

t  thickness  (m) 

U  gas  velocity  vector  (ms-1) 

v  gas  velocity  (ms-1) 

V  voltage  (V) 

Ujy  volume  occupied  by  Nation  in  catalyst  layer  (m3 ) 

vs  volume  occupied  by  solid  catalyst  particles  in  cata¬ 

lyst  layer  (m3) 

vt  total  volume  of  catalyst  layer  (m3 ) 
vv  volume  of  void  space  in  catalyst  layer  (m3 ) 

Y  mass  fraction 

X  mole  fraction 

x,  y,  z  Cartesian  coordinates 

Greek  Symbols 

a  charge  transfer  coefficient  for  the  reaction 

X  Bruggemann  coefficient 

£gdl  porosity  of  gas  diffusion  layer 

£N,cat  volume  fraction  of  Nation  in  catalyst  layer 
£s,cat  volume  fraction  of  solid  catalyst  particles  in  catalyst 
layer 

£v,cat  volume  fraction  of  void  space  in  catalyst  layer 
0  electric  potential  (V) 

cp  phase  potential  (V) 

A  gas  channel  width  ratio  {Id Ip) 

A  membrane  water  content  (kmol  H20(kmol  S03  -  )-1 ) 

p  density  of  gas  (kg  m-3 ) 

a  electronic  conductivity  of  solid  catalyst  particle 

(£2-1  m-1) 

r  ionic  conductivity  (£2-1  m-1 ) 

r  fluid  stress  (N  m-2) 

£  source  term  in  electronic  conduction  equation 

Superscripts 
eff  effective 

Subscripts 
a  anode 

b  bottom 

c  cathode/gas  channel 

Cat  catalyst  layer 

GDI  gas  diffusion  layer 

in  inlet 


m  membrane 

N  Nation 

P  carbon  plate 

out  outlet 

ref  reference 

S  solid  catalyst  particles 

t  top 

1-6  Cell  1 -Cell  6 


However,  in  general,  the  voltage  produced  by  a  single-cell  fuel 
cell  only  ranges  from  0.65  to  1.0  V.  For  practical  applications,  a  num¬ 
ber  of  cells  must  be  assembled  in  series  to  construct  a  fuel  cell  stack. 
The  number  of  unit  cells  required  in  a  stack  depends  on  the  desired 
voltage  or  the  net  power  capacity,  which  can  vary  from  less  than 
100  W  for  portable  power  generators  to  around  300  kW  for  heavy 
vehicles  [7].  Recently,  significant  progress  in  the  development  of 
PEMFC  stack  technology  has  been  achieved  by  an  increasing  num¬ 
ber  of  experimental  [8-11]  and  theoretical  [12-16]  studies. 

When  designing  a  fuel  cell  stack,  one  needs  be  able  to  mea¬ 
sure  individual  cell  voltages.  In  an  ideal  situation,  all  unit  cells  are 
identical,  and  there  are  no  performance  differences  between  the 
cells.  In  reality,  cell  performances  differ,  sometimes  drastically,  due 
to  unequal  reactant  supply.  In  practice,  non-evenly  distributed  cell 
voltages  among  the  unit  cells  may  cause  a  great  reduction  in  the 
power  output  of  the  stack.  To  achieve  an  even  voltage  distribution, 
the  stack  should  be  operated  on  the  premise  that  the  reactant  gases 
are  distributed  uniformly  to  all  unit  cells  of  the  stack.  Since  both 
the  fuel  and  the  oxidant  gases  are  fed  to  the  anode  and  the  cathode 
sides,  respectively,  of  all  unit  cells  of  the  stack  through  a  manifold 
and  the  flow  channel  network,  the  design  of  the  manifold  and  the 
flow  channel  network  is  regarded  as  one  of  the  critical  issues  related 
to  stack  performance.  Karimi  et  al.  [17]  found  that  the  stack  perfor¬ 
mance  is  considerably  affected  by  the  size  and  the  number  of  flow 
channels  on  bipolar  plates.  To  achieve  improved  stack  performance, 
the  stack  manifold,  flow  channels  on  the  bipolar  plates,  and  the 
stack  operating  conditions  should  be  designed  carefully.  Mohamed 
and  Jenkins  [12]  presented  an  approach  based  on  the  genetic  algo¬ 
rithm  technique  to  find  the  optimal  configuration  of  the  stack  in 
terms  of  the  number  of  cells  and  cell  surface  area.  Pan  et  al.  [18] 
carried  out  an  optimal  design  of  complex  manifold  geometries  for 
uniform  flow  distribution  between  microchannels  by  means  of  the 
flow  network  method.  Tonomura  et  al.  [19]  performed  a  CFD-based 
optimal  design  of  a  manifold  in  plate-fin  microdevices. 

In  the  present  study,  numerical  methods  are  employed  to  inves¬ 
tigate  the  dependence  of  the  performance  of  a  six-cell  PEMFC  stack 
on  the  geometric  design  of  the  gas  channels  in  the  individual  cells. 
Adjustment  for  the  channel  width  ratios  (zl)  or  the  heights  of  gas 
channels  (h)  of  the  cells  has  been  attempted  to  obtain  evenly  dis¬ 
tributed  cell  voltages.  In  the  meantime,  a  six-cell  fuel  cell  stack 
(190  W)  has  been  assembled  for  experiments,  and  the  polarization 
curves  of  the  stack  have  been  adopted  to  verify  the  validity  of  the 
numerical  solutions. 

2.  Theoretical  model 

The  physical  model  of  a  stack  consisting  of  six  unit  cells,  each  of 
which  is  equipped  with  parallel  straight  channels  on  both  the  anode 
and  the  cathode  carbon  plates,  is  shown  in  Fig.  1.  There  are  sev¬ 
eral  advantages  associated  with  the  straight  channels.  For  example, 
compared  to  the  serpentine  channels,  where  there  may  be  eddies 
near  bends  in  the  channel,  straight  channels  provide  fewer  places 
for  water  to  accumulate  in  the  channels.  Another  advantage  is  that 
there  is  a  smaller  pressure  drop  with  fluids  flowing  in  the  channels 
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Computation 
domain 


(a)  Six-cell  stack  with  parallel  straight  channels 


(b)  Computation  domain  containing  six  cells 


Fig.  1.  Schematic  of  a  six-cell  proton  exchange  membrane  fuel  cell  stack,  (a)  Six-cell 
stack  with  parallel  straight  channels,  (b)  Computation  domain  containing  six  cells. 


since  there  are  no  corners.  Furthermore,  plates  with  straight  par¬ 
allel  channel  flow  fields  may  be  easier  to  manufacture  than  those 
with  more  complex  shapes.  Fig.  1(a)  illustrates  a  schematic  of  the 
stack  with  the  cells  numbered  from  1  to  6,  and  Fig.  1(b)  conveys  the 
computation  domain  of  the  solution  module,  which  covers  a  half 
channel  area  for  the  six  cells  and  the  inlet/outlet  manifolds.  For  the 
pattern  of  straight  parallel  channels  considered,  if  the  number  of 
gas  channels  is  sufficiently  large,  a  periodic  transport  phenomenon 
is  expected  to  develop  from  channel  to  channel  for  all  the  cells. 
Thus,  one  only  needs  to  deal  with  the  computation  domain  indi¬ 
cated  by  the  dashed  lines,  and  the  right  and  the  left  faces  ( y-z 
planes  at  x  =  0  and  x  =  lP)  of  the  computation  domain  are  the  sym¬ 
metric  planes.  Flere,  the  origin  of  the  Cartesian  coordinate  is  put 
at  the  bottom  left  corner  of  the  computation  domain  where  the  y- 
direction  represents  the  streamwise  direction  of  the  channels.  The 
reactant  gases,  hydrogen  and  air,  are  delivered  into  the  cells  through 
the  inlet  manifolds,  and  the  product  gases  exit  through  the  outlet 
manifolds  after  the  reaction.  For  each  cell,  the  gas  channels  on  the 
anode  side  and  those  on  the  cathode  side  are  placed  in  a  staggered 
arrangement.  The  width  of  the  solution  domain  is  denoted  with  lP, 
and  the  half  channel  width  is  lc:  therefore,  the  channel  has  width 
ratio  (A)  defined  as  lc/lp.  The  channel  width  and  channel  height 
(h)  are  regarded  as  two  major  geometric  parameters  to  elevate  the 
performance  of  the  stack.  The  values  of  the  fixed  parameters  are 
listed  in  Table  1.  The  existing  information  for  the  values  of  these 
fixed  parameters  is  scattered.  However,  the  parameters  are  chosen 
based  on  a  review  of  the  existing  information  [1-6].  Note  that  in 
this  study,  the  operating  temperature  and  pressure  are  assigned  to 


Table  1 

Fixed  parameters. 


Parameter 

Symbol 

Value 

Operating  temperature,  K 

T 

343 

Operating  pressure,  Pa 

P 

101,325 

Operating  stack  voltage,  V 

V 

3.5-5.75 

Gas  channel  length,  m 

k 

5.0  x  10-2 

Width  of  a  solution  module,  m 

Ip 

1.0  x  10“3 

Thickness  of  catalyst  layer,  m 

teat 

3.0  x  10“5 

Thickness  of  carbon  plate,  m 

tp 

3.0  x  10-3 

Thickness  of  gas  diffusion  layer,  m 

tGDL 

3.0  x  10“4 

Thickness  of  membrane,  m 

bn 

1.76  x  10“4 

Bruggemann  coefficient  of  gas  diffusion  layer 

Xgdl 

1.5 

Bruggemann  coefficient  of  catalyst  layer 

XCat 

1.5 

GDL  permeability,  m2 

1<GDL 

1.76  x  10-11 

Catalyst  layer  permeability,  m2 

h-Cat 

1.76  x  10-11 

Anode  inlet  hydrogen  velocity,  m  s-1 

Va,in 

0.1667 

Cathode  inlet  air  velocity,  m  s-1 

Vc,in 

0.5 

Anode  inlet  mass  fraction  of  H2 

nw 

0.2 

Anode  inlet  mass  fraction  of  H20 

^H20,i 

0.8 

Cathode  inlet  mass  fraction  of  02 

^o2)i 

0.174 

Cathode  inlet  mass  fraction  of  N2 

^N2)i 

0.609 

Cathode  inlet  mass  fraction  of  H20 

^H20,i 

0.217 

Anode  concentration  parameter 

ra 

1.0 

Cathode  concentration  parameter 

rc 

1.0 

Anodic  charge  transfer  coefficients  for  anode 

da, a 

0.5 

reaction 

Cathodic  charge  transfer  coefficients  for  anode 

°ta,c 

0.5 

reaction 

Aanodic  charge  transfer  coefficients  for  cathode 

Ole, a 

1.5 

reaction 

Cathodic  charge  transfer  coefficients  for  cathode 

Olc,c 

1.5 

reaction 

Reference  transfer  current  density  at  anode,  Am-3 

ja,ref 

4.5  x  109 

Reference  transfer  current  density  at  cathode, 

Jc.ref 

250 

Am-3 

Volume  fraction  of  Nafion  in  catalyst  layer 

SN,Cat 

0.3 

Volume  fraction  of  void  space  in  catalyst  layer 

£v,Cat 

0.112 

Volume  fraction  of  solid  catalyst  particles  in 

Ss.Cat 

0.588 

catalyst  layer 

Porosity  of  gas  diffusion  layer 

£gdl 

0.5 

Effective  ionic  conductivity  of  catalyst  layer, 

pejf 

Cat 

4.2 

£2_1  m_1 

Ionic  conductivity  of  Nafion,  m_1 

r  N.Cat 

25.56 

Effective  mass  diffusivity  of  oxygen  in  catalyst 

Deff 

LJ02,Cat 

1.9546  x  10“3 

layer,  cm2  s_1 

Effective  mass  diffusivity  of  hydrogen  in  catalyst 

Deff 

H2,Cat 

0.985  x  10“3 

layer,  cm2  s_1 

Effective  mass  diffusivity  of  oxygen  in  GDL,  cm2  s_1 

Deff 

o2,gdl 

Deff 

h2,gdl 

1.845  x  10-2 

Effective  mass  diffusivity  of  hydrogen  in  GDL, 

0.930  x  10-2 

cm2  s_1 

Electronic  conductivity  of  gas  diffusion  layer, 

or  GDL 

150 

£2_1  ITT1 

Electronic  conductivity  of  carbon  plate,  m_1 

Op 

4000 

Electronic  conductivity  of  catalyst  layer,  £2_1  m_1 

&S,Cat 

1414 

be  343  K  and  101,325  Pa,  respectively.  The  value  of  the  operating 
voltage  of  the  fuel  cell  stack  is  varied  between  3.5  and  5.75  V  in  the 
computation. 

The  theoretical  model  is  developed  under  the  following  simpli¬ 
fications  and  assumptions: 


1.  Compression  pressure  is  applied  by  fastening  the  bolts  that 
clamp  the  cells  of  the  stack  tightly,  and,  therefore,  at  the  inter¬ 
faces  between  the  layers  of  any  cells,  there  is  no  appreciable 
contact  resistance  against  the  conduction  of  heat,  electrons,  or 
ions. 

2.  For  all  the  cells,  the  gas  flows  in  the  gas  channels,  the  porous 
layers  are  laminar  and  compressible,  and  the  effects  of  buoyancy 
on  the  flow  motion  are  ignored. 

3.  The  gas  diffusion  layers  and  the  catalyst  layers  are  homogeneous, 
isotropic  porous  media. 
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4.  Inlet  and  outlet  conditions  and  thermodynamic  and  electro¬ 
chemical  properties  of  the  gases  and  the  solid  materials  of  the 
stack  components  are  assumed  to  be  constant. 

5.  Water  exists  in  the  entire  fuel  cell  only  in  the  vapor  phase,  and 
no  phase  change  process  is  taken  into  consideration. 

6.  The  operating  temperature  is  uniform  and  constant  in  the  stack. 

7.  Because  of  the  high  ohmic  resistance  of  the  membrane,  the  inter¬ 
face  between  the  catalyst  layer  and  the  membrane  is  assumed  to 
be  insulated  against  electrons. 

8.  The  fuel  cell  stack  is  operated  under  steady  conditions. 

The  present  three-dimensional  model  involves  electrochemical 
reactions;  mass,  momentum,  and  species  transport  phenomena; 
ionic  conduction;  electronic  conduction;  and  Nation  loading  effects. 
The  modeling  proposed  by  Cheng  et  al.  [5]  that  was  applied  for  a 
single-cell  system  is  extended  to  the  analysis  of  the  six-cell  stack. 
The  governing  equations  of  the  theoretical  model  are  described 
briefly  in  the  following. 

2.1  Electrochemical  reactions  (catalyst  layers) 

Note  that  for  each  cell,  H2  and  02  are  consumed  by  the  reactions 
in  the  catalyst  layers  of  the  anode  and  cathode  sides,  respectively, 
and  H20  is  the  product  of  the  reaction  at  the  cathode  sides.  The 
source  terms  for  different  species  are 


r  J  a 

%2=“2 F 

(la) 

00 

o 

II 

*n\n 

(lb) 

NJ 

o 

II 

1 

° 

(lc) 

where  ja  and  jc  are  the  transfer  current  densities  at  the  anode 
and  the  cathode,  respectively,  which  are  calculated  based  on  the 
Butler-Volmer  condition: 


Momentum  conservation: 

V  •  (sflUU)  -  -sVP  +  V  ■  (et)  +  {ELfJL  (4) 

K 

Species  equation: 

V(spUYi)  =  V-(DfvYi)  +  Si  (5) 

where  £  and  k  are  the  porosity  and  the  permeability  of  the  porous 
layers,  respectively;  r  represents  the  fluid  stress  tensor.  Fuel  gas 
at  the  anode  contains  H2  and  H20,  and  oxidant  gas;  the  cathode 
contains  02,  N2,  and  H20.  Y;  denotes  the  species  mass  fraction  of 
species  i.  Note  that  for  pure  fluid  flows  in  the  gas  channels,  the 
values  of  s  and  k  are  s  =  1  and  k  =  oo. 

In  Eq.  (5),  the  term  S;  standing  for  the  source  term  of  species  i 
is  calculated  by  Eq.  (1),  and  the  calculation  for  the  effective  mass 
diffusivity  of  species  i  (Dff )  is  based  on  Bruggemann’s  equation, 
that  is, 

Df  =  (e)xDj  (6) 

where  Dz  is  the  mass  diffusivity  of  species  i,  and  x  is  the  Bruggemann 
coefficient  of  the  porous  layer. 

2.3.  Ionic  conduction  (catalyst  layers  and  proton  exchange 
membranes) 

Ionic  conduction  occurs  in  the  catalyst  layers  and  the  proton 
exchange  membrane.  The  ionic  conduction  equation  for  determi¬ 
nation  of  the  phase  potential  (<p)  in  the  catalyst  layers  is 

-V-(r%tV(pCat)  =  $  (7) 

where  (pCat  is  the  phase  potential  and  r^t  is  the  effective  ionic  con¬ 
ductivity  of  the  catalyst  layers.  The  effective  ionic  conductivity  r^t 
is  evaluated  by  using 


Ja  —Jajef 


jc  —Jcjef 


Yh7 


*H2 ,  ref 


“p{  "}] 


y<h 


Y02,ref 


exp{^}-exp{^,}] 


RT 

—dc,cF 


(2a) 

(2b) 


where  ja >re/and  jcref  are  the  reference  transfer  current  densities,  and 
ra  and  rc  the  concentration  parameters,  at  the  anode  and  the  cath¬ 
ode,  respectively.  In  addition,  aa,a  and  aa,c  are  the  anodic  and  the 
cathodic  charge  transfer  coefficients  for  anode  reaction,  and  ac,a , 
and  aCtC  the  anodic  and  the  cathodic  charge  transfer  coefficients 
for  cathode  reaction.  Furthermore,  Yo2jef  and  YH2jef  are  the  refer¬ 
ence  concentrations  at  which  the  exchange  current  densities  were 
obtained.  All  the  reference  values  of  these  coefficients  are  given  in 
Table  1. 


2.2.  Mass,  momentum,  and  species  transport  (gas  channels,  GDLs, 
and  catalyst  layers) 


The  flow  channels  and  the  gas  diffusion  layers  in  the  cells  are 
employed  to  deliver  the  reactant  gases  uniformly  into  the  catalyst 
layers  for  reactions.  The  governing  equations  to  predict  the  mass, 
momentum,  and  species  transport  phenomena  in  the  gas  channels, 
GDLs,  and  catalyst  layers  can  be  expressed,  in  general  form,  in  accor¬ 
dance  with  the  conservation  laws  in  mass,  momentum,  and  species. 
In  the  porous  medium  layers,  the  equations  are  derived  based  on 
the  non-Darcy  law. 

Mass  conservation: 


V  •  (epU)  =  0 


nl  =  ^%tr"'Cat  (8) 

In  addition,  for  the  catalyst  layers,  the  source  term  of  Eq.  (7),  £,  is 
equal  to  the  transfer  current  densities  -ja  and  jc  for  the  anode  and 
the  cathode,  respectively,  which  are  calculated  with  Eqs.  (2a)  and 
(2b). 

On  the  other  hand,  the  ionic  conduction  equation  for  the  proton 
exchange  membrane  is  expressed  as 


v-(rmv<?m)  =  o 


(9) 


where  (pm  is  the  phase  potential  and  rm  is  the  ionic  conductivity  of 
the  membrane,  which  can  be  determined  by  the  relation  proposed 
by  Springer  et  al.  [20]  as 


rm{T)  =  exp 


,268  (a® 


)  (0.005239A.  —  0.00326) 


(10) 


where  the  water  content  of  the  membrane  (A)  is  dependent  on 
a  dimensionless  water  activity  in  the  membrane  (a)  according  to 
experimental  data.  The  activity  is  defined  as 


a  = 


Xh2qP 

Psat 


(11) 


where  the  saturation  pressure  of  water  vapor  can  be  computed 
from 


log wPsat  =  -2.1794 +  0.02953(1 -273.15) -9.1837 

X  10“5(T  -  273. 15)2  +  1 .4454  x  10“7(T  -  273. 15)3 


(3) 


(12) 
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with  T  in  K.  The  correlation  between  A  and  a  is  derived  as 

f  0.043  +  17.81a  -  39.85a2  +  36a3  0  <  a  <  1 
A  =  <  (13) 

1^14  +  1.4(0-1)  1  <  a  <  3 

Based  on  the  solutions  for  the  phase  potential  (<p),  the  distri¬ 
bution  of  local  ionic  current  density  (i)  in  the  catalyst  layers  and 
the  proton  exchange  membrane  may  be  further  determined  with 

i  =  -rv<p. 

2.4.  Electronic  conduction  (carbon  plates,  GDLs,  and  catalyst 
layers) 

The  carbon  plate,  the  GDL,  and  the  catalyst  layers  all  serve  as 
conductors  for  electric  current.  In  the  carbon  plate  and  the  GDL, 
there  is  no  electrochemical  reaction  taking  place.  The  electrochem¬ 
ical  reactions  only  take  place  in  the  catalyst  layers.  Therefore,  in 
these  three  components,  the  electronic  conduction  equation  may 
be  expressed  in  general  form  as 

V  •  {aeffV(p)  =  f  (14) 

where  <p  is  the  electronic  potential  and  oeH  denotes  the  effective 
electronic  conductivity.  For  the  carbon  plates,  the  effective  elec¬ 
tronic  conductivity  ae^  is  treated  as  a  constant  (crP).  For  the  porous 
gas  diffusion  and  catalyst  layers,  the  magnitude  of  the  effective 
electronic  conductivity  is  determined  in  terms  of  the  porosity  in 
accordance  with  the  Bruggemann’s  equation  as 

^l  =  ^-£cdl)XgdloCDl  (15a) 

acat  =  £s%as,cat  (15b) 

It  is  noticed  that  the  effective  electronic  conductivity  of  the  cat¬ 
alyst  layer  (or^t)  increases  with  the  volume  fraction  of  solid  catalyst 
particles  in  the  catalyst  layer  ( £s,cat )•  Also,  note  that,  for  the  catalyst 
layers,  the  source  term  f  is  identical  to  that  of  Eq.  (7),  whereas  for 
the  carbon  plates  or  the  GDL,  the  source  term,  f ,  is  zero  for  without 
an  electrochemical  reaction. 

2.5.  Nafion  loading  (catalyst  layers) 

When  Nafion  is  introduced  into  the  catalyst  layer,  it  reduces  the 
electronic  conductivity  of  the  catalyst  layer.  However,  since  it  is 
a  good  conductor  for  the  hydrogen  protons,  Nafion  loading  can 
elevate  the  ionic  conductivity  of  the  catalyst  layers.  In  addition, 
the  introduction  of  Nafion  into  the  catalyst  layers  simultaneously 
changes  the  porosity  and  the  gas  permeability  of  the  porous  layers 
so  that  the  transport  behavior  of  the  gas  species  is  affected. 

As  Nafion  is  loaded,  the  total  volume  of  the  catalyst  layer  (vt) 
may  be  divided  into  three  major  portions,  namely,  the  volume  occu¬ 
pied  by  the  solid  catalyst  particles  (Us),  the  volume  occupied  by 
Nafion  (vN),  and  the  volume  of  the  void  space  (tty).  The  volume 
fractions  of  the  three  portions,  appearing  in  Eqs.  (8)  and  (15),  are 
denoted  by  £s,cat  =  vS/Vt,  £N,cat  =  vN/Vt,  and  £v,Cat  =  l-£cfcat“£jv,cat. 
respectively. 

2.6.  Boundary  conditions 

1.  At  the  inlet,  the  inlet  velocities  for  the  anode  gas  (hydrogen)  and 
the  cathode  gas  (air)  in  the  gas  channel,  Va  in  and  VCt  in,  are  fixed  at 
0.166  and  0.5  ms-1,  respectively.  The  stoichiometries  of  hydro¬ 
gen  and  air  corresponding  to  the  inlet  velocities  are  1.25  and 
3.85  at  V=4.2  V,  respectively.  The  inlet  mass  fractions  for  all  the 
gas  species  are  listed  in  Table  1.  Note  that  both  the  hydrogen 
and  the  air  enter  the  fuel  cell  at  100%  relative  humidity.  The  flow 
boundary  condition  is  fixed  independent  of  the  cell  performance, 


implying  that  this  is  a  constant  flow  operating  mode  compared 
to  the  constant  stoichiometry  mode. 

2.  At  the  outlet,  the  exit  pressure  is  prescribed  to  be  101,325  Pa. 

3.  The  width  of  a  module  (Ip)  is  set  to  be  1  mm.  Therefore,  for  the 
faces  of  the  solution  domain  at  x  =  0  and  x=l  mm,  symmetric 
boundary  conditions  for  U,  Y,  0,  and  cp  are  applied. 

4.  The  solid  walls  of  the  gas  channels  on  the  carbon  plates  are 
impermeable,  and  hence,  on  the  solid  walls,  3yI/3n  =  0  and  U  =  0. 

5.  The  ionic  conduction  takes  place  only  inside  the  catalyst  lay¬ 
ers  and  proton  exchange  membranes,  and  the  protons  are  not 
allowed  to  transfer  into  other  components;  therefore,  on  the 
faces  of  the  catalyst  layers  and  the  membranes  the  normal  gra¬ 
dient  of  phase  potential  is  assigned  to  be  zero. 

6.  Since  both  the  gas  channels  and  the  membranes  are  insulators 
against  electrons,  at  the  boundaries  of  the  gas  channels  and  the 
membranes,  the  normal  gradient  of  the  electric  potential  (0) 
vanishes. 


This  set  of  conservation  equations,  along  with  the  property 
equations,  is  solved  by  adopting  a  finite-volume  scheme  on  struc¬ 
tured  grids  within  the  framework  of  the  commercial  computational 
fluid  dynamics  code  (CFD-ACE+).  The  number  of  control  volumes 
varies  from  502,850  to  710,084  to  check  the  grid-independence  of 
the  numerical  predictions,  and  the  comparison  in  the  average  cur¬ 
rent  density  has  been  made  based  on  the  results  obtained  from 
different  grids.  It  is  found  that  for  the  case  at  V=4.2V,  the  rela¬ 
tive  error  is  within  1.2%.  Therefore,  typically  the  grid  containing 
502,850  control  volumes  is  used  to  save  computation  time  without 
loss  of  accuracy.  In  computation,  the  iteration  for  the  solutions  for 
U,  Y,  0,  and  cp  progresses  until  the  following  convergence  criterion 
is  satisfied: 


10“ 


(16) 


where  ^  represents  each  of  the  variables,  U,  Y,  0,  and  cp\  l  is  the 
iteration  number;  and  i,j  and  k  represent  the  spatial  indices  for  the 
x-,  y-  and  z-directions,  respectively. 

The  average  current  density  (1)  of  fuel  cells  is  calculated  by  inte¬ 
grating  a  local  electric  current  density  distribution  on  the  carbon 
plate  surface  at  z  =  0  of  cell  1  divided  by  the  carbon  plate  surface 
area  as  I  =  1  /A  is  ■  dA,  where  is  is  the  local  electric  current  den¬ 
sity  on  the  carbon  plate  surface,  and  A  is  the  carbon  plate  surface 
area. 

It  is  important  to  note  that  condensation  of  water  vapor  may  be 
critical  to  the  performance  of  the  six-cell  stack.  Nevertheless,  for  the 
present  case,  experiments  show  that  the  water  condensation  at  the 
cathode  side  is  not  severe,  compared  to  the  problem  of  non-evenly 
distributed  channel  flow  rate  at  V=  4.2  V,  even  in  a  long-term  oper¬ 
ation.  This  probably  may  be  attributed  to  the  over-stoichiometric 
operating  condition.  The  stoichiometry  of  air  is  set  to  be  3.85  so 
that  the  air  stream  may  definitely  help  remove  the  water  conden¬ 
sate.  Meanwhile,  the  straight  channels  might  also  significantly  aid 
the  water  condensate  removal.  However,  at  a  lower  stoichiometry, 
or  for  other  types  of  channels  where  the  operating  voltage  is  lower, 
the  effects  of  water  condensation  vapor  will  become  severe.  In  that 
case,  assumption  (5)  may  not  be  valid.  On  the  other  hand,  the  effects 
of  the  temperature  distribution  may  play  an  important  role  in  the 
condensation  and  the  vaporization  of  water.  Therefore,  the  temper¬ 
ature  field  in  the  fuel  cell  stack  should  also  be  included  because  the 
condensation  of  water  vapor  is  severe.  The  conditions  for  the  liquid 
water  to  exist  were  investigated  extensively  by  Wang  and  cowork¬ 
ers  [21,22].  Assumptions  (5)  and  (6)  can  be  partially  justified  with 
the  observation  from  the  experimental  and  computational  works 
on  the  two-phase  flow  in  PEMFC  of  Refs.  [21,22]. 
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(a)  Experimental  apparatus 


(b)  Six-cell  stack  module  (190-W) 

Fig.  2.  Experimental  apparatus  and  six-cell  stack  module,  (a)  Experimental  appara¬ 
tus  (b)  Six-cell  stack  module  (190  W). 


3.  Experiments  with  six-cell  PEMFC  stack 

A  six-cell  stack  (190W)  with  cross-sectional  area  of 
22  cm  x  22  cm  and  active  area  of  14  cm  x  14  cm  has  been  installed 
and  tested  to  partly  confirm  the  present  approach’s  numerical 
predictions.  The  purpose  of  the  experiments  is  simply  to  test  the 
validity  of  the  commercial  code.  The  operating  and  geometric 
conditions  in  the  simulation  are  consistent  with  the  experimental 
module.  Seven  composite  graphite  plates  of  3-mm  thickness  are 
used  as  the  bipolar  plates.  The  parallel  straight  gas  channels 
machined  on  the  bipolar  plates  are  1  mm  wide  and  1  mm  deep. 
Six  5-layer  membrane  exchange  assemblies  (MEA)  provided  by 
Dupont  Inc.  are  placed  at  the  center  between  two  bipolar  plates  in 
each  cell.  The  membranes  used  in  this  study  are  Nation®  112,  and 
the  gas  diffusion  layers  are  made  of  carbon  fiber  papers.  Before 
entering  the  channels,  the  fuel  and  oxidant  gases  are  humidified  to 
maintain  100%  relative  humidity  and  343  K  temperature. 

In  addition,  to  obtain  the  polarization  curves  with  the  fuel  cell, 
a  performance  test  system  has  also  been  set  up.  In  the  experi¬ 
ments,  the  volumetric  flow  rates  of  the  hydrogen  and  the  air  are 
determined  such  that  the  inlet  gas  velocities  agree  with  the  values 
assigned  in  the  numerical  simulation,  that  is,  Vain  =  0.1667 ms-1 
and  l/bin  =  0.5ms_1.  The  fuel  cell  temperature  is  maintained  at 
343  K.  Photographs  of  the  experimental  apparatus  and  the  six-cell 
stack  are  shown  in  Fig.  2. 

A  comparison  between  the  experiments  and  the  numerical  pre¬ 
diction  is  required  for  validation  of  the  model.  Actually,  a  number  of 
experiments  have  been  performed  for  comparison.  Herein,  only  a 
typical  comparison  is  shown  to  avoid  duplication  and  for  the  sake  of 


Fig.  3.  Polarization  curves:  comparison  between  predictions  and  experiments. 


brevity.  The  comparison  between  the  numerical  predictions  and  the 
experimental  data  for  the  polarization  curves  of  the  six-cell  stack 
module  is  provided  in  Fig.  3.  Note  that  the  reference  current  den¬ 
sities,  transfer  coefficients  and  reference  concentrations  have  been 
specified  so  as  to  fit  the  experimental  data.  It  has  been  observed 
that  the  simulation  model  can  be  tuned  to  closely  agree  with  exper¬ 
imental  data.  In  this  figure,  it  is  found  that  at  V=4.2  V  (average  cell 
voltage  is  0.7  V),  the  numerical  and  the  experimental  data  are  both 
approximately  2600  Am-2.  In  this  study,  the  numerical  simulation 
is  typically  performed  for  the  cases  at  4.2  V. 

4.  Results  and  discussion 

A  parametric  study  of  the  effects  of  the  channel  width  ratios 
(yl)  and  the  channel  heights  (ft)  has  been  carried  out,  and  a  better 
combination  of  the  size  of  the  gas  channels  for  the  individual  cells 
can  be  determined  for  evenly  distributed  cell  voltages.  In  this  study, 
seven  test  cases  with  different  distributions  of  channel  width  ratios 
and  channel  heights  are  investigated,  which  are  listed  in  Table  2.  In 
this  table,  A\  to  Aq  represent  the  channel  width  ratios,  and  h\  to 
ft6  the  channel  heights,  for  cells  1  to  6,  respectively.  Among  these 
cases,  Case  (1)  is  the  base-line  case  used  in  the  experimental  stack 
system,  in  which  all  the  sizes  of  the  cells  are  the  equal.  In  Cases  (2) 
to  (7),  the  sizes  of  the  channels  differ  in  each  of  the  different  cells. 
The  relative  performance  of  these  cases  will  be  compared  with  Case 
(1). 

Fig.  4  shows  the  voltage  and  the  current  density  distributions 
at  different  x-z  planes  for  Case  (1).  In  the  case  shown  in  Fig.  4(a), 


Table  2 

Test  cases  with  different  combinations  of  channel  width  ratios  and  channel  heights. 
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the  zero  voltage  is  set  at  the  bottom  surface  of  the  stack  (Vb  =  0), 
and  the  voltage  at  the  top  surface  (Vt)  is  assigned  to  be  4.2  V.  It  is 
seen  that  in  the  stack,  the  voltage  distribution  in  the  computation 
domain  is  decreased  from  4.2  to  0  V.  On  the  other  hand,  the  distri¬ 
butions  of  the  current  density  at  each  of  these  same  planes  exhibit 
a  periodic  pattern  from  cell  to  cell.  However,  along  the  channels, 
it  is  observed  that  the  strength  of  the  current  density  is  gradually 
reduced.  This  may  be  attributed  to  the  consumption  in  the  hydro¬ 
gen  and  oxygen  flowing  in  the  channels  that  causes  a  decrease  in 
the  local  concentrations  of  the  gases  downstream. 

Mass  fraction  distributions  of  oxygen  at  x  =  0.75  mm  (central 
cross  section  on  y-z  plane  on  the  cathode  channels  side)  for  Case 
(1 )  at  different  operating  voltages  are  displayed  in  Fig.  5.  In  this  fig¬ 
ure,  cases  at  V=  3.6, 4.2,  and  4.8  V  are  plotted.  In  general,  the  highest 
oxygen  concentrations  are  observed  in  the  inlet  manifold  area  and 
the  lowest  in  the  outlet  manifold  area.  It  is  found  that  at  V=4.8  V, 
higher-concentrations  of  oxygen  prevail  in  the  channels  of  all  cells 
and  that  the  distribution  of  oxygen  concentration  is  also  uniform. 
However,  as  the  operating  voltage  is  reduced  to  3.6  V,  the  con¬ 
centration  of  oxygen  in  the  cells  decreases  remarkably,  and  in  the 
downstream  zone,  the  supply  of  oxygen  is  obviously  insufficient. 
Furthermore,  it  is  observed  that  the  distribution  of  mass  fraction  of 
oxygen  in  the  cells  is  uniform  at  a  higher  voltage  (e.g.,  4.8  V);  how¬ 
ever,  at  lower  voltages,  the  equality  of  oxygen  concentration  in  the 
cells  is  no  longer  visible.  Note  that  the  more  uniform  distribution 
at  higher  cell  voltages  results  from  the  constant  flow  rate  operat¬ 
ing  mode.  Since  there  is  much  less  consumption  of  reactant  gases 
at  lower  current  densities,  the  distribution  of  concentration  and 
thus  performance  is  more  uniform.  On  the  other  hand,  because  the 
reactant  flow  is  kept  constant,  at  a  lower  voltage  (associated  with 
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(a)  Voltage  distribution 
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Fig.  4.  Voltage  and  current  density  distributions  at  different  x-z  planes  for  Case  ( 1 ). 
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Fig.  5.  Mass  fraction  distributions  of  oxygen  at  x  =  0.75  mm  (central  cross  section  on 
the  cathode  channels  side),  for  Case  (1)  at  different  operating  voltages. 


higher  current  densities)  it  corresponds  to  a  lower  stoichiometry 
situation.  Therefore,  at  lower  voltage  operations,  the  oxygen  con¬ 
centration  is  not  uniform.  At  V=3.6  V,  the  concentration  of  oxygen 
in  Cells  3  and  4  is  appreciably  lower  than  that  in  other  cells,  and 
Cells  1  and  6  possess  highest  concentration  of  oxygen.  Also  note 
that  the  cell  voltage  results  from  the  gas  concentration  in  the  cells. 
Thus,  one  may  expect  that  Cells  1  and  6  possess  higher  cell  voltages 
than  Cells  3  and  4.  Fig.  6  demonstrates  the  predicted  cell  voltages  for 
individual  cells  of  the  stack  at  various  operating  voltages  for  Case 
(1).  It  is  clearly  found  that  the  individual  cell  voltages  at  V=4.8V 
are  relatively  uniform,  whereas  at  V=  3.6  V  the  individual  cell  volt¬ 
ages  are  obviously  different  for  different  cells.  As  expected,  Cells 
1  and  6  possess  higher  cell  voltages,  and  Cells  3  and  4  have  lower 
voltages.  In  accordance  with  the  obtained  results,  it  is  recognized 
that  an  even  voltage  distribution  results  from  a  uniform  reactant 
gas  distribution  over  the  flow  fields  of  the  cells.  Therefore,  in  this 
study,  the  channel  width  ratios  or  the  heights  of  gas  channels  are 
adjusted  for  individual  cells  to  yield  a  uniform  distribution  of  the 
reactant  gases.  It  is  noted  that  the  sizes  of  the  channels  are  the  same 
for  both  the  cathode  and  the  anode  bipolar  plates.  Therefore,  only 
the  cathode  channels  sizes  are  mentioned  in  the  remainder  of  the 
paper. 

4.1.  Adjustment  in  channel  width  ratios 

As  stated  in  the  preceding  section,  Case  (1)  has  with  some  dis¬ 
advantages,  and  there  is  still  room  for  significant  improvement  in 
the  uniformity  of  cell  voltages,  especially  at  lower  operating  volt¬ 
age.  Therefore,  various  new  design  models  are  further  tested  in  this 
study  by  numerical  experiments  to  determine  a  better  design. 
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Cell  number 

Fig.  6.  Cell  voltages  for  individual  cells  of  the  stack  at  various  operating  voltages,  for 
Case  (1). 


Figs.  7  and  8  show  the  mass  fraction  distributions  of  oxygen  at 
the  same  cross  section  for  Case  (2)  and  Case  (3),  respectively,  at 
different  operating  voltages.  In  Case  (2),  to  increase  the  oxygen  con¬ 
centration  in  Cells  3  and  4  and  decrease  the  oxygen  concentration 
in  Cells  1  and  6,  one  elevates  the  channel  width  ratio  of  Cells  3  and 
4  to  be  0.55  and  also  reduces  the  width  ratio  of  Cells  1  and  6  to 
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Fig.  7.  Mass  fraction  distributions  of  oxygen  at  x  =  0.75  mm  (central  cross  section  on 
the  cathode  channels  side),  for  Case  (2)  at  different  operating  voltages. 
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Fig.  8.  Mass  fraction  distributions  of  oxygen  at  x  =  0.75  mm  (central  cross  section  on 
the  cathode  channels  side),  for  Case  (3)  at  different  operating  voltages. 


be  0.45.  Likewise,  the  flows  in  Cells  3  and  4  are  indeed  accelerated, 
and  those  in  Cells  1  and  6  are  retarded.  However,  contrary  to  the 
observations  in  Case  (1),  in  Case  (2)  the  concentration  of  oxygen  in 
Cells  3  and  4  becomes  higher  than  that  in  other  cells,  and  Cells  1 
and  6  possess  a  lower  concentration  of  oxygen.  The  distribution  of 
oxygen  concentration  is  still  not  uniform  over  the  flow  fields  of  the 
cells.  In  Case  (3),  the  channel  width  ratio  of  Cells  3  and  4  is  further 
elevated  to  0.6,  and  that  of  Cells  1  and  6  is  reduced  to  be  0.4.  In  Fig.  8, 
it  is  seen  that  in  this  case  the  uniformity  becomes  much  poorer  and 
hence  not  acceptable. 

According  to  the  results  with  Cases  (1),  (2),  and  (3),  it  is  rea¬ 
sonable  to  expect  a  better  design  that  is  a  compromise  between 
Cases  (1 )  and  (3).  Therefore,  after  a  series  of  computations  for  var¬ 
ious  size  combinations,  the  channel  width  ratios  are  adjusted  to 
be  A\  =0.495,  A2  =  0.500,  A3  =  0.508,  A4  =  0.508,  yl5  =  0.502,  and 
As  =  0.498,  while  the  heights  of  channels  for  all  cells  are  maintained 
at  1  mm.  This  case  is  referred  to  as  Case  (4)  as  listed  in  Table  2.  The 
relative  performance  of  Case  (4)  is  shown  in  Fig.  9,  as  compared 
with  Cases  (1)  and  (2).  The  distributions  of  average  mass  fraction 
of  oxygen  cross  at  y  =  28.5  mm  (at  middle  point  of  channel  length). 
A  comparison  of  the  uniformity  of  concentration  is  made  for  these 
three  cases  at  a  lower  operating  voltage,  3.6  V.  In  this  figure,  it  is 
obvious  that  the  design  of  Case  (4)  leads  to  a  more  uniform  concen¬ 
tration  distribution  in  the  cells.  Detailed  information  regarding  the 
concentration  distribution  in  the  six  cells  for  Case  (4)  is  presented 
in  Fig.  10.  In  this  figure,  distributions  of  mass  fraction  of  oxygen 
at  the  interfaces  between  membranes  and  cathode  catalyst  layers 
in  the  six  cells  for  Case  (4)  are  plotted.  It  is  observed  that  the  dis¬ 
tributions  in  the  six  cells  are  nearly  identical.  This  means  that  the 
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Fig.  9.  Distributions  of  average  mass  fraction  of  oxygen  on  the  cross  section  at 
y  =  28.5  mm  (at  middle  point  of  channel  length)  in  the  cells,  for  Cases  (1 ),  (2),  and  (4) 
at  V  =  3.6  V. 
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Fig.  10.  Distribution  of  mass  fraction  of  oxygen  at  the  interfaces  between  mem¬ 
branes  and  cathode  catalyst  layers  in  the  six  cells  for  Case  (4)  at  V=3.6  V. 
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Fig.  11.  Mass  fraction  distributions  of  oxygen  atx  =  0.75  mm  (central  cross  section  on 
the  cathode  channels  side),  for  Case  (5)  at  different  operating  voltages,  (a)  V=4.8  V 
(b)  V=4.2V  (c)  V  =  3.6V. 


designed  manifold  is  able  to  supply  the  reactant  gases  to  the  six 
cells  evenly. 

It  is  important  to  mention  that  an  even  distribution  of  cell  volt¬ 
ages  is  also  found  with  Case  (4)  when  used  together  with  the 
uniform  concentration  distribution.  As  a  result,  the  power  density 
with  Case  (4)  is  increased  by  16.5%  at  V=4.2V  (1  =  3030  Am-2),  as 
compared  to  the  value  with  Case  (1 ). 


4.2.  Adjustment  in  channel  heights 

It  has  been  recognized  that  channel  height  could  be  one  of  the 
influential  geometric  parameters  on  the  distribution  of  the  reactant 
gases  in  the  cells  [5,6].  Next,  the  channel  heights  in  the  six  cells  are 
adjusted  to  improve  the  uniformity  of  the  reactant  gas  distribution. 

Cross-sectional  mass  fraction  distributions  of  oxygen  at 
x  =  0.75  mm  for  Cases  (5)  and  (6)  at  different  operating  voltages  are 
displayed  in  Figs.  11  and  12,  respectively.  In  Case  (5),  for  increas¬ 
ing  the  oxygen  concentration  in  Cells  3  and  4  and  decreasing  the 
oxygen  concentration  in  Cells  1  and  6,  the  channel  heights  of  Cells 
3  and  4  are  increased  to  1.05  mm,  and  the  channel  heights  of  Cells 
1  and  6  are  reduced  to  0.95  mm.  In  these  two  cases,  the  channel 
width  ratio  for  all  the  cells  is  fixed  at  0.5.  Fig.  11  shows  that,  in  Case 
(5),  the  flows  in  Cells  3  and  4  are  over-accelerated,  while  those  in 
Cells  1  and  6  are  over-retarded.  As  a  result,  the  concentration  of 
oxygen  in  Cells  3  and  4  becomes  higher  than  that  in  other  cells,  and 
Cells  1  and  6  now  possess  lower  concentrations  of  oxygen.  This  is 
definitely  not  the  desired  situation. 


358 


C.-H.  Cheng,  H.-H.  Lin  /  Journal  of  Power  Sources  194  (2009)  349-359 


z  =  26.028  mm 


z  =  0  mm 

y  =  0  mm  (a)  V=  4.8  Volts 


(c)  V=  3.6  Volts 


Present  cross 
section 
(x=0.75  mm) 

1 

x 

y  =  67  mm 


Fig.  12.  Mass  fraction  distributions  of  oxygen  at  x  =  0.75  mm  (central  cross  section 
on  the  cathode  channels  side),  for  Case  (6)  at  different  operating  voltages. 


Furthermore,  considering  Case  (6)  in  which  the  channel  heights 
of  Cells  3  and  4  are  further  elevated  to  1.2  mm,  and  those  of  Cells 
1  and  6  are  reduced  to  be  0.8  mm,  one  observes  in  Fig.  12  that  the 
uniformity  of  gas  concentration  becomes  worse,  especially  at  low 
operating  voltages.  For  example,  at  V=3.6V,  the  oxygen  concen¬ 
trations  in  Cells  1  and  6  are  rather  low,  indicating  an  insufficient 
oxygen  supply  for  the  two  cells. 

Similarly,  the  results  presented  in  Figs.  5, 11  and  12  reveal  that  a 
better  combination  for  the  channel  heights  for  different  cells  should 
be  a  compromise  between  Cases  (1)  and  (5).  Therefore,  a  series  of 
computations  for  various  size  combinations  has  been  carried  out  in 
order  to  determine  the  best  combination  for  the  channel  heights  for 
different  cells.  The  suggested  design  is  listed  in  Table  2  and  num¬ 
bered  as  Case  (7).  In  Case  (7),  the  channel  heights  are  adjusted  to 
hi  =  0.98  mm,  h2  =  1.0  mm,  h3  =  1.02  mm,  h4  =  1.02  mm,  h5  =  1.0  mm, 
and  h6  =  0.98mm,  while  the  channel  width  ratio  is  fixed  at  0.5. 
Fig.  13  conveys  the  distributions  of  the  average  mass  fraction  of 
oxygen  on  the  cross  section  at  y  =  28.5  mm  (at  the  mid-point  of  the 
channel  length)  in  the  cells  for  Cases  (1),  (5),  and  (7)  at  V=3.6V.  In 
this  figure,  it  is  noticed  that  the  design  of  Case  (7)  results  in  a  more 
uniform  concentration  distribution  in  the  cells.  As  compared  with 
Case  (1),  the  improvement  in  gas  concentration  distribution  in  the 
cells  with  Case  (7)  results  in  an  increase  in  the  power  density  by 
13%  at  V= 4.2  V  (1  =  2940  A  m-2 ). 

It  may  be  worth  mentioning  that  the  optimization  of  the  channel 
width  ratio  (and  height)  distribution  is  dependent  on  the  mani¬ 
fold  design,  the  physical  properties,  and  the  operation  parameters. 
While  the  size  distribution  has  been  found  to  be  optimal  for  this  cell 
design,  this  is  not  a  universal  distribution  for  all  manifolds  or  stack 
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Fig.  13.  Distributions  of  average  mass  fraction  of  oxygen  on  the  cross  section  at 
y  =  28.5  mm  (at  middle  point  of  channel  length)  in  the  cells,  for  Cases  (1 ),  (5),  and  (7) 
at  V=  3.6  V. 
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Table  3 

Electric  current  density  at  various  operating  voltages  for  the  test  cases. 


designs.  Also  note  that  the  pressure  drop  associated  with  the  chan¬ 
nel  design  is  an  important  issue.  Flowever,  here  in  the  present  study 
while  the  sizes  of  some  channels  are  increased,  the  sizes  of  some 
other  channels  are  decreased  so  that  the  average  channel  width 
and  height  are  kept  at  A  =  0.5  and  h  =  1.0.  In  addition,  in  the  paral¬ 
lel  straight  channels,  the  pressure  drop  with  fluids  flowing  in  the 
channels  is  much  less  compared  to  other  types  of  channels.  There¬ 
fore,  the  variation  of  pressure  drop  is  not  addressed  here.  Table  3 
shows  the  numerical  data  for  the  electric  current  density  at  various 
operating  voltages  for  the  test  cases.  As  expected,  Case  (4)  exhibits 
a  relatively  higher  current  density  at  both  V=  3.6  and  4.2  V  as  com¬ 
pared  to  other  test  cases,  whereas  Case  (2)  leads  to  a  lower  current 
density  due  to  the  non-uniform  distribution  of  the  reactants  in  the 
cells. 

5.  Concluding  remarks 

In  this  study,  a  numerical  simulation  of  the  transport  phe¬ 
nomena  in  a  six-cell  PEM  fuel  cell  stack  has  been  performed  to 
investigate  the  possibility  of  obtaining  evenly  distributed  cell  volt¬ 
ages  by  adjusting  the  channel  width  ratios  (yl)  or  the  heights  of 
gas  channels  (h)  of  the  cells.  Meanwhile,  a  six-cell  fuel  cell  stack 
(190  W)  has  been  assembled  for  experiments,  and  the  polarization 
curves  of  the  stack  have  been  adopted  to  verify  the  validity  of  the 
numerical  solutions. 

In  the  base-line  case,  Case  ( 1 ),  it  is  found  that  at  V=  4.8  V,  higher- 
concentrations  of  oxygen  prevail  in  the  channels  of  all  the  cells, 
and  the  distribution  of  oxygen  concentration  is  also  uniform.  How¬ 
ever,  as  the  operating  voltage  is  reduced  to  3.6,  or  even  2.8  V,  the 
concentration  of  oxygen  in  the  cells  decreases  remarkably,  and  in 
the  downstream  zone,  the  supply  of  oxygen  is  obviously  insuffi¬ 
cient.  Furthermore,  it  is  observed  that  the  concentration  of  oxygen 
in  Cells  3  and  4  is  appreciably  lower  and  higher  in  Cells  1  and  6. 
Note  that  under  the  assumption  of  no  dry  out  in  the  cells,  an  even 
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voltage  distribution  can  be  obtained  if  the  distribution  of  reactant 
gas  is  uniform  among  the  cells.  Therefore,  in  this  study,  the  channel 
width  ratios  and  the  heights  of  gas  channels  are  adjusted  to  yield  a 
uniform  distribution  of  the  reactant  gases  in  the  cells. 

Results  show  that  an  even  distribution  in  cell  voltages  can  be 
obtained  simply  by  adjusting  the  channel  width  ratios  for  individ¬ 
ual  cells.  Based  on  the  numerical  simulation,  it  is  found  that  the 
design  of  Case  (4)  leads  to  a  more  uniform  cell  voltage  distribution 
in  the  cells.  In  Case  (4),  the  channel  width  ratios  are  A\  =0.495, 
A2  =  0.500,  A3  =  0.508,  A4  =  0.508,  A5  =  0.502,  and  A6  =  0.498, 
while  the  heights  of  channels  for  all  cells  are  maintained  at  1  mm. 
Based  on  the  design  of  Case  (4),  the  power  density  is  increased  by 
16.5%  at  V=4.2  V,  as  compared  to  the  value  with  Case  (1). 

On  the  other  hand,  an  even  distribution  in  cell  voltages  can  also 
be  obtained  by  adjusting  the  channel  heights  for  individual  cells. 
A  better  combination  of  the  channel  heights  for  individual  cells, 
which  is  labeled  as  Case  (7),  has  been  found,  and  the  combination  is 
h\  =  0.98  mm,  h2  =  1.0  mm,  h3  =  1.02  mm,  h4  =  1.02  mm,  ft5  =  1.0  mm, 
and  h6  =  0.98  mm,  while  the  channel  width  ratio  is  fixed  at  0.5.  As 
compared  with  Case  (1 ),  the  improvement  in  gas  concentration  dis¬ 
tribution  in  the  cells  with  Case  (7)  results  in  an  increase  in  the  power 
density  by  13%  at  V= 4.2  V. 
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